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Abstract 

In this paper we discuss how to generate inductive invariants for safety verifi- 
cation of hybrid systems. A hybrid symbohc-numeric method is presented to 
compute inequahty inductive invariants of the given systems. A numerical 
invariant of the given system can be obtained by solving a parameterized 
polynomial optimization problem via sum-of-squares (SOS) relaxation. And 
a method based on Gauss-Newton refinement and rational vector recovery 
is deployed to obtain the invariants with rational coefficients, which exactly 
satisfy the conditions of invariants. Several examples are given to illustrate 
our algorithm. 

Keywords: semidefinite programming, sum-of-squares relaxation, safety 
verification, invariant generation 



1. Introduction 

Complex physical systems are systems in which the techniques of sensing, 
control, communication and coordination are involved and interacted with 
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each other. Among complex physical systems, many of them are safety crit- 
ical systems, such as airplanes, railway, and automotive applications. Due 
to the complexity, ensuring correct functioning of these systems, e.g., spatial 
separation, especially collision avoidance of aircrafts during the entire flights, 
is among the most challenging and most important problems in computer sci- 
ence, mathematics and engineering. 

As a common mathematical model for complex physical systems, hybrid 
systems j^, [sl are dynamical systems that are governed by interacting dis- 
crete and continuous dynamics [11, @, 0]. Continuous dynamics is specified 
by differential equations, which is possibly subject to domain restrictions or 
algebraic relations resulting from physical circumstances or the interaction of 
continuous dynamics with discrete control. For discrete transitions, the hy- 
brid system changes state instantaneously and possibly discontinuously, for 
example, the instantaneous change of control variables like the acceleration 
(e.g., the changing of a by setting a := —h with braking force > 0). 

The verification of hybrid systems is an important problem that has been 
studied extensively both by the control theory, and the formal verification 
community for over a decade. Among the most important verification ques- 
tions for hybrid systems are those of safety, i.e., deciding whether a given 
property if) holds in all the reachable states, and the dual of safety, i.e., reach- 
ability, deciding if there exists a trajectory starting from the initial set that 
reaches a state satisfying the given property ip. In principle, safety verifica- 
tion or reachability analysis aims to show that all trajectories of the hybrid 
systems starting from the initial set cannot enter some unsafe regions in the 
state space. 

Safety verification or reachability analysis of hybrid systems presents a 
more difficult challenge, primarily due to the infinite number of possible states 
in continuous state space. Some well-established techniques have been pro- 
posed. In [2I, [ll|, quantifier elimination was used to calculate exact reachable 
sets for linear systems with certain eigenstructures and semialgebraic initial 
sets. Tiwari |28| generalized this method to handle linear systems with al- 
most arbitrary eigenstructures. In level set methods, ellipsoidal 
techniques and flow-pipe approximations have been presented for computing 
approximate reachable sets of hybrid systems. 

Recently, some methods 20|,l2l|,|25|,|28|] based on invariant generation have 



been proposed for safety verification of hybrid systems. An invariant [2J] of 
a hybrid system is a property that holds in all the reachable states of the 
system, in other words, it is an over- approximation of all the reachable states 
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of the system. Invariants are useful facts about the dynamics of a given 
system, and are widely used in numerous approaches to analyze and verify 
systems. For example, if the invariants lie inside the safe regions, or their 
intersection with the unsafe regions is empty, then safety of hybrid systems 
is verified. 

The problem of generating invariants of an arbitrary form is known to 
be computationally hard, intractable even for the simplest classes. The 
usual technique for generating invariants is to produce an inductive invari- 
ant, i.e., an assertion that holds at the initial states of the system, and is 
preserved by all discrete and continuous state changes. There has been a 
considerable volume of work towards invariant generation for hybrid sys- 
tems using techniques in convex optimization, semi-algebraic system solving 
0, H [isl, [3, [H M, H 0, H [3, liol- However, some of these techniques 



are only applicable to linear systems, some are subject to numerical errors 
and some suffer from high complexity. In virtue of the efficiency of nu- 
merical computation and the error-free property of symbolic computation, a 
hybrid symbolic-numeric method via SOS relaxation and exact certificate is 
presented in 3l| to construct inequality invariants for continuous dynamic 
systems given by nonlinear vector fields. 

In this work, we study how to generate inequality invariants for safety ver- 
ification of nonlinear hybrid systems. We present a hybrid symbolic-numeric 
method, based on sum-of-squares (SOS) relaxation via semidefinite program- 
ming (SDP) and exact SOS representation recovery, to generate inequality 
invariants of hybrid systems, which guarantee that all the reachable states 
never enter the given unsafe regions. The idea is as follows: (1) Given a safe 
property, we predeterminate the templates of the invariants, and construct 
a semidefinite programming (SDP) system to solve the corresponding para- 
metric polynomial optimization problem. (2) An exact invariant is obtained 
by recovering the exact SOS representation from the approximate solution 
of the associated SDP system. In the recovery step, Gauss-Newton iteration 
is deployed to refine the approximate solution from SDP solver. Then safety 
property of the hybrid systems can be easily verified, by the exact SOS rep- 
resentations of the conditions of the invariants. More details will be shown 
in Section [31 

Unlike the numerical approaches, our method can yield exact invariants, 
which can overcome the unsoundness in the verification of hybrid systems 



caused by numerical errors [19[. In comparison with some symbolic ap- 



proaches of invariant generation based on qualifier elimination technique. 
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our approach is more efficient and practical, because parametric polynomial 
optimization problem, based on SOS relaxation method, can be solved in 
polynomial time theoretically. 

The rest of the paper is organized as follows. In Section |2l we introduce 
some notions about hybrid systems and invariants. Section |3] is devoted to 
illustrating a symbolic-numeric approach to generate invariants for safety 
verification of hybrid systems. In Section m we present some examples on 
invariant generation for safety verification of hybrid systems. Section |5] con- 
cludes the paper and discusses some future work. 

2. Invariants 

To model hybrid systems, we recall the definition of hybrid automata j^, 

isj- 

Definition 1 (Hybrid system). A hybrid system H : {V, L,T,Q,T),'^ , io) 
consists of the following components: 

• V = a set of real-valued system variables. A state is an 
interpretation ofV, assigning to each Xi & V a real value. An assertion 
is a first-order formula over V. A state s satisfies an assertion ip, 
written as s \= ip, if (f holds on the state s. We will also write tpi \= ip2 
for two assertions fi,(p2 to denote that ip2 is true at least in all the 
states in which Lpi is true; 

• L, a finite set of locations; 

• T, a set of (discrete) transitions. Each transition r : {£,£', g^-, p^-) G 
T consists of a prelocation i E L, a postlocation i' G L, the guard 
condition g-r over V , and an assertion pr over V VJ V representing 
the next-state relation, where V = denotes the next-state 
variables. Note that the transition r can take place only if g^. holds; 

• G, an assertion specifying the initial condition; 

• V, a map that maps each location i E L to a differential rule (also 
known as a vector field or a flow field^ ^(^)> of the form Xi = fe,i{V) 
for each Xi G V , written briefly as = f£(x). The differential rule at a 
location specifies how the system variables evolve in that location; 
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• a map that maps each location i & L to a location condition (loca- 
tion invariant) "^{i), an assertion overV; 

• £0 & L, the initial condition. We assume that the initial condition 
satisfies the location invariant at the initial location, that is, 6 |= \E'(^o)- 

By a state of a hybrid system H : (V, L, T, O, T), ,£q), we mean the tuple 
{£, x) G L X M" where n is the number of program variables in H. 

Definition 2 (Computation). fi3 / A computation of a hybrid system H is 
an infinite sequence of states 

< /qi >, < /i, Xi >,■■■,< /j, Xj >, < /j+i, Xj+i >, • • • 

such that 

• [Initiation] /q = £0 and xq |= 6; 

Furthermore, for each consecutive pair < U^^i >, < /i+i,Xj+i >, one of 
the two consecution conditions holds: 

• [Discrete Consecution] There exists a transition r : {£,£', g^-, p^-) 
such that li = /j+i = £' and (xj,Xj+i) |= p^(xj, Xj+i) if g^- holds, or 

• [Continuous Consecution] /j = /j+i = £, and there exists a time 
interval 6 > and a smooth (continuous and differentiahle to all or- 
ders) function f : [0,(5] — )■ s.t. f evolves from Xj to Xj+i according 
to the differential rule 'D{£) at location £, while satisfying the location 
invariant '^{£). Formally, 

- /(O) = x„ f{6) = x,+i and Wt e [0, S],f{t) \= ^{£), 

- VtG [0,5),(/(t),/(t)) h^^W- 

A state {£, x) is a reachable state of a hybrid system H if it appears in a 
computation of H. 

Figure [1] is a graphical representation of a hybrid system with two loca- 
tions £i,£2- A state of this hybrid system is denoted by {£, x) G {£i,£2} xW^, 
and the initial state set is ^1 x 0. During a continuous flow, the discrete 
location £i is maintained and the continuous state variables x evolve accord- 
ing to the differential equations x = /^^(x), with x satisfying the location 
invariant '^{£i). At the state (^j,x), if the guard condition g{£i,£j) is met, 
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Figure 1: An example of hybrid system H 

the system may undergo a transition to location ij, and x will take the new 
value x', which is determined by the reset map p{£i, £j). 

Given a hybrid system with an initial set and a prespecified safe (or 
unsafe) region, the system is safe if starting from any state in the initial 
set, this system would never evolve to the given unsafe region or the system 
would always stay inside the safe region. More specifically, consider the 
hybrid system H shown in Figure [1] and let Xu C M" be an unsafe region. 
The system H is said to be safe if all trajectories of the system starting from 
any state in (^i,Xo) G £i x 6 can not reach X^, or any state in X^ is not 
reachable. 

In this work, we will apply the invariant generation method to verify 
safety of hybrid systems. The following definitions of invariants of hybrid 
systems come from (25| . 

Definition 3 (Invariant). An invariant of a hybrid system at location i is 
an assertion X such that for any reachable state {£, x) of the hybrid system, 

X h^:. 

An invariant of a hybrid system is an assertion that holds in all the reach- 
able states of the system. 

The problem to generate invariants with arbitrary form is known to be 
computationally hard, intractable even for the simplest classes. The usual 
technique for generating invariants is to compute inductive invariants, defined 
as follows. 

Definition 4 (Inductive invariant). An inductive assertion map X of a hy- 
brid system H : {V, L,T,Q,'D,'$ , £q) is a map that associates with each lo- 
cation i E L an assertion X{i) that holds initially and is preserved by all 
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discrete transitions and continuous flows o/H. More formally, an inductive 
assertion map satisfies the following requirements: 



(i) [Initial] ^I{io). 

(ii) [Discrete Consecution] For each discrete transition r : {£, i' , Qr, Pr) , 

starting from a state satisfying and taking r leads to a state sat- 
isfying X{i') . Formally, 

I{e)AgrApr |=X(f) 

where represents the assertion with the current state vari- 

ables Xi, . . . ,Xn replaced by the next state variables x[, . . . , x'^, respec- 
tively. 

(iii) [Continuous Consecution] For every location £ G L and states (£, xi), 
(£, X2) such that'x.2 evolves from 'x.^ according to the differential rule T>[pj 
at I, if xi \= X{1) then Xa |= . 

Our definition of inductive invariants is sliglitly modified from tliat of 



Definition 4 in (25[, and ttie only change made is taking the guard conditions 
into account. 

For a hybrid system, a formula X(x) is called a differential invariant at 
location i if X(x) satisfies conditions (i) and (iii), that is, X(x) holds initially 
and is preserved by the continuous flow at a single location. There are sev- 
eral literature to compute differential invariants. [l3| presented an approach 
based on the computable algebraic-geometry theory to generate differential 



invariants. 



20| computed differential invariants using a verification logic for 



hybrid systems. [3l| suggested a hybrid symbolic-numeric method to com- 



pute inequality differential invariants. 

Remark 1. Clearly, inductive invariants are over- approximation of the reach- 
able sets of hybrid systems, since an inductive invariant is true for all the 
reachable states of the system. 



3. Safety Verification of Hybrid Systems 

The aim of this section is to translate the problem of safety verification 
of hybrid systems into that of generating invariants, which can be trans- 
formed further into polynomial optimization problem with parameters. We 



7 



will present a hybrid symbolic-numeric method, based on SOS relaxation, to 
solve this polynomial optimization problem, and obtain the invariants, which 
can guarantee the safety property of hybrid systems. 

3.1. Invariants and Safety Verification 

In this paper, we are interested in hybrid systems in which the relations 
are given by (real) polynomials over the system variables. Then we define 

Definition 5 (Polynomial Hybrid System). A polynomial hybrid system is 
a hybrid system: H : (V, L, T, 6, I^, £o); where 

• for each transition r : {£,£', gr, Pr) G T, the guard condition gr (resp. 
the reset relation p^) is a conjunction of polynomial inequalities over V 
(resp. V U V); also, the initial condition and the location invari- 
ant "^{i), for each £ E L, are conjunctions of polynomial inequalities 
over V; 

• each rule 'D{t) is of the form Xi = /f,j(x) for each Xi G V , where 

/,,,(x) eM[x]. 

We are interested in finding invariants of the form (/^^(x) > at loca- 
tion £ G L. Below is an alternative expression of Definition |H 

Theorem 1. Let H : (y, L, T, 6,'D,\l/,^o) be a hybrid system. Suppose 
for each location i G L, there exists a function (pe{'x) satisfying the following 
conditions: 

(i) (x)>o, 

(ii) V5^(x) > A g{i, i') A p{i, i') |= (fieix.') > 0, for any transition {i, f , g, p) 

going out from i, 

(iii) (/^^(x) > A \E'(£) 1= (/3£(x) > 0, here (peipi) denotes the Lie-derivative of 
Lpe along the vector field V{i), i.e., (piijs.) = Yl^^^i fff ■ fe,i{^), 

then (feipi.) > is an invariant of the hybrid system H at location i. 

Proof. The proof follows directly from Definition |H □ 

Remarked that if the functions ^pei'x.) at all locations are identical to </3(x), 
then V3(x) is an inductive invariant of the given hybrid system, as described 
in the following theorem. 
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Theorem 2. Let H be a hybrid system. Suppose there exists a function v?(x) 
satisfying the following conditions: 

(i) e h vi^) > 0, 

(ii) v?(x) > A g{i,i') A p{ij') 1= ifiix.') > 0, for any transition {£,£', g,p) 
going out from i, 

(iii) (/9(x) > A 1= v?(x) > 0, 

then ip{'x) > is an inductive invariant of the system H. 

In the sequel, for brevity, we shall use ifei'x.) to denote both the invari- 
ant fei'x.) > and ipeijx.). 

The following theorem shows that invariants can be applied to verify the 
safety property of hybrid systems. 

Theorem 3. Let H be a hybrid system, and Xu{i) be the unsafe region at 
location i. Suppose there exists functions V9^(x), for £ ^ L, that satisfy the 
conditions (i-iii) in TheoremUl and moreover, 

(iv) Xu{i) h M^) < 0, V£ G L, 

then the safety of the system H is guaranteed. 

Proof. Clearly, (feipi.) > is an invariant of hybrid system H at location i. 
Then the condition (iv) implies that all reachable sets at location £ lie outside 
the unsafe region Xu{£), yielding the safety of the system. □ 

Similarly, inductive invariants can be applied to verify safety of hybrid 
systems. 

Theorem 4. Let H be a hybrid system, and Xu{i) be the unsafe region at 
location £. Suppose there exists a function (/?(x) that satisfies the conditions 
(i-iii) in Theorem\^ and moreover, 

(iv) X^{£) h ¥'(x) < 0, V£ e L, 

then the safety of the system is guaranteed. 

Proof. The proof is similar to that of Theorem [3l □ 

Remark 2. Functions (p/fx) and V5(x) in Theorems\^ and^ are also known 
as barrier certificates in \2n l. 
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3.2. Sum of Squares Relaxation 

According to Theorems [3] and HI to verify the safety of hybrid system H, 
it suffices to compute real polynomials (/^^(x) or <^(x). 

In the following, we only discuss how to find the invariant (/^^(x) at each 
location £ e L. The problem of computing the inductive invariant (/'(x) can 
be handled similarly. 

Our idea of computing </?£(x) or (/^(x), based on Sum-of-Squares (SOS) 
relaxation and rational vector recovery, is as follows. 

Step 1: Predetermine a template of polynomial invariants with the given 
degree and convert the problem of computing polynomial invariants 
to the associated parametric polynomial optimization problem. SOS 
relaxation method is then applied to obtain a polynomial invariant 
with fioating point coefficients. 

Step 2: Apply Gauss- Newton refinement and rational vector recovery on 
the approximate polynomial invariant to get polynomials with rational 
coefficients, which exactly satisfy the conditions of invariants of the 
given hybrid system. 

The problem of computing the invariant (/^^(x) at each location £ G L, that 
satisfy the conditions in Theorem [3] can be transformed into the following 
problem 



Let us first predetermine a template of polynomial invariants with the 
given degree d, that is, we assume 



where x" = ■ ■ ■ x"" and Cq, G M are parameters, with a G Zi>q and 
^iLi (^i ^ d. One can apply quantifier elimination methods to solve the 
corresponding parametric semi-algebraic systems, and for the given tem- 
plate, quantifier elimination methods can yield the sufficient and necessary 



find v^^(x) G M[x],V£ G L 



s.t. eh¥'£oW>0, 

Vii.-^) > A g{l, f ) A f ) h ^i' (x') > 0, > 

(/p,(x) > A ^(^) h > 0, 
Xu{i) 1= M^) < 0. 



(1) 




(2) 



a 
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conditions for the existence of invariants. Several Maple packages, such as 



RAGLib|l6| and DISCOVERER [32|, are available to solve this problem. 



However, quantifier elimination method based on the cylindrical algebraic 
decomposition (CAD) are of high complexity. Instead, we will explore the 
SOS relaxation techniques based on semidefinite programming (SDP) solving 
to obtain polynomial invariants. 
In the sequel, we suppose that 

q p 

e = {x G : /\ ^Kx) > 0}, X^{i) = {x e M" : /\ 0,,(x) > 0}, 
1=1 j=i 



^{i) = {x e M" : /\ ^,,fc(x) > 0}, g{i,i') = {x G M" : /\ .(x) > 0}, 

k=l i=l 
t 

p(£,f)(x,x') = {x' e : /\ p«',.(x,x') > 0}, 



u=l 



where G L, and 6'i(x), V^£,fc(x), gu',i{?^) and p«/_„(x,x') are poly- 

nomials. 

Clearly, a sufficient condition for r(x) G ]R[x] with degree 2e to be positive 
semidefinite is that there exists an SOS of r(x): 

r(x) = ^ri2(x), with ri(x) G M[x], (3) 

i 

or, equivalently, r(x) can be represented as 

r(x) = m(x)^ ■ W ■ m(x), 

where ly is a real symmetric and positive semidefinite matrix, and m(x) is 
a vector of terms in ]R[x] with degree < e. 

When a polynomial r(x) can be written as an SOS in M[x], we simply call 
r(x) an SOS. Denote by T,n,2e the set of all SOSes of degree < 2e in variables 
X\ , . . . , I.e. , 

S„,2e = {r(x) G M[x] : r(x) is an SOS, deg(r(x)) < 2e}. 

Based on the SOS relaxation, the constraints in ([T]) can be replaced by 
stronger ones. For instance, to find a polynomial v^£o(x) satisfying 

e |=<^,„(x) >o 
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it suffices to find (x) such that 

(/?£o(x) = ao(x) + ^ a,(x)0i(x), 
1=1 

where o"o,a; G M[x] are SOSes. Therefore, the problem of computing poly- 
nomials V5£(x) is transformed into the following SOS program: 

find <^^(x) e M[x],V£ G L 

s.t. V5£o W = c^o(x) + Ya=i crii^)Oi{x), 

V9£/(x') = A«/,o(x) + Xli=i A«',i(x)5(«/,i(x) I 
+ ELi7^£',«(x)p«',„(x,x') + w(x)(/3^(x), I 
ipe{x) = 0£,o(x) + ^fc^i 0£,fc(x)V'£_fc(x) + i/£(x)v?£(x) + ££,1 
-(^^(x) = /i£,o(x) + ^J^i /x^j(x)Oj(x) + e£,2, 

where (Ti(x), A«/,i(x), 7«/,t.(x), ?7«/(x), 0^,fc(x), //^(x), /i^j(x) G E„,2e and e£,i, 
Q,2 £ The decision variables are the coefficients of all polynomials 

appearing in (jl]), such as (/^^(x), (T;(x), A£^/_i(x). 

Since the coefficients of ipi{x.),r]£gi{x.) and ^'^(x) are unknown, some non- 
linear terms that are products of these coefficients, occur in the second and 
third constraints of (jlj). The SOS relaxation will then lead to a non-convex 
bilinear matrix inequalities (BMI) problem. To avoid BMI problem, we adopt 
stronger conditions to compute the invariants of hybrid systems. 

Theorem 5. Under the assumptions in TheoremUl suppose for each i E L, 
^Piipi.) satisfies the following conditions: 

(i) eh¥'^o(x)>o, 

(ii') g{i,i') A p{i,i') \= '^I'i'x.') > 0, for any transition {£,£', g,p) going out 
from i, 

(iii') ^(i) 1= v3,(x) > 0, 

then ifiipi.) > is an invariant of the hybrid system H at location i. In 
addition, if Lp£{x.) satisfies 

(iv) Xu{i) hv'Kx) <0,V£gL, 

then the safety of the system is guaranteed. 
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Proof. Since the conditions (ii') and (iii') are stronger than the conditions 
(ii) and (iii) in Theorem [T] respectively, (fg is an invariant at location L 
According to Theorem [3l the condition (iv) can guarantee the safety of this 
system. □ 

A similar conclusion can be attained for inductive invariants, as stated in 
the following 

Theorem 6. Under the assumptions in Theorem suppose there exists a 
polynomial y^i'x) satisfying the following conditions: 

(i) e 1= (^(x) > 0, 

(ii') g{i,i') A p{i,i') 1= V^(x') > 0, for any transition {£,£', g,p) going out 
from i, 

(iii') h V5(x) > 0,V^ e L, 

then V2(x) > is an inductive invariant of the hybrid system H. In addition, 
if (p{x.) satisfies 

(iv) Xu{i) hv'lx) <0,V^eL, 

then the safety of the system is guaranteed. 

Having Theorem [5l the program can be modified into the following 
problem: 

find (^^(x) G M[x],V£ G L 

s.t. (fioi^) = o-o(x) + XlLi M^Wii^)^ 

(/9^/(x') = A«/,o(x) + J2i=i >^u'A^)9u',i{^) + Yll=i 7«',«(x)p«',«(x,x'), 

ipeix.) = 0£,o(x) + Ylk=i 0^,fc(x)i/'£,fc(x) + e£,i, 

-V9£(x) = ^£,o(x) + /^^j (x)Oj(x) + e£,2, 

(5) 

where cr;(x), A«/,i(x), 7a/,„(x), ^^^^(x), //£j(x) G S„,2e and e^,!, ee,2 G K+- The 
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program is equivalent to the following SDP problem: 



inf Tr(M, W, V, P, Q) 

s.t. ^,„(x) = mo(x)^ ■ MM ■ mo(x) + ^Li ^^x)^ ■ ■ mKx)^;(x), 



+ Ell w«',(x)^ ■ ■ W«, ,(x)(7«,,(x) 

+ Yfu=i^u'A^V ■ "1^'"''"' ■ v«/_„(x)pa/,„(x,x'), 
¥.,(x) = P.,o(x)^ ■ P[^'°] ■ p,,o(x) 

+ ELi PiA^f ■ P^'''^ ■ P^,fc(x)V^,,fc(x) + e,,i, 
(/..(x) = -q,,o(x)^ ■ gt^'O] ■ q,,o(x) 

- Ej=i q^,j(x)^ ■ Q'^'-''' ■ q£j(x)C£j(x) - e^,2, 

(6) 

where all the matrices M^'l, W^^^ y^^^ Qt^'-'l are symmetric and pos- 

itive semidefinite, and the function Tr(M, W, V, P, Q) denotes the sum of 
traces of all these matrices, which acts as a dummy objective function com- 
monly used in SDP for optimization problem with no objective function. 



Many Matlab packages of SDP solvers, such as SOSTOOLS [2^, YALMIP 



14| . and SeDuMi 26|], are available to solve the problem ([6]) efficiently. 



3. 3. Exact Certificate of Sum of Squares Decomposition 

Since the SDP solvers in Matlab is running in fixed precision, the tech- 
niques in Section 13.21 will yield numerical solutions to the associated SDP 
problem ([6]), where the numerical polynomial 'Pei'x.) and numerical positive 
semidefinite matrices M^''\ . . . , Q^^'^^ satisfy the constraints in ([6]) approxi- 
mately, for instance, 

<^,„(x) ^ mo(x)^ ■ MM ■ mo(x) + ^Li ^^x)^ ■ ■ mKx)^;(x), ^ 0. 

(7) 

However, due to round-off errors, (/^^(x) > may not necessarily be an 
invariant of the given hybrid system at location i, because the constraints 
in (ini) may not hold exactly, for example, (I7j) means that feoi'x.) may not 
be positive semidefinite exactly within the initial set O. Therefore in the 
next step, from the numerical polynomials ^pti'x.) and the numerical positive 
semidefinite matrices M''', . . . , Q^^'^\ we will recover polynomials ^^(x) with 
rational coefficients, which satisfy ([6]) exactly. 

As described in [9|, finding a polynomial with rational coefficients can 
be translated into the problem of rational vector recovery. In Section 13. 2^ 
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a numerical vector v| denoting the coefficient (column) vector of (pe{'x.) is 
obtained by solving the SDP system, i.e., (/^^(x) = vl"'" ■ Ti{x), where T^^x.) is 
the column vector of all terms in (/^^(x). To obtain a rational vector near 
to v^, we can employ the simultaneous Diophantine approximation algorithm 



12l |. once the bound of the common denominator of is given. 

The recovery of the matrices Aft'^ , . . . , Q^^'^^ into rational positive semidefi- 
nite matrices is split into two steps. We ffist recover the matrices M^, . . . , Q'^'-'' 
for l<l<q,...,l<j<p and then recover MM,...,g[^'°l To illustrate 
the idea, we only discuss how to recover M^'^ for 1 < I < q and the matri- 
ces W^^^ V^^^ pt^''^], Q^^'^^ can be recovered similarly. 

Given the numerical positive semidefinite matrices M''', 1 < I < q in 
we can find the nearby rational positive semidefinite matrices M''' by use of 
the rational vector recovery technique. In practice, all the M^'l are numerical 
diagonal matrices, in other words, the off-diagonal entries are very tiny and 
the diagonal entries are nonnegative. Therefore, by setting the small entries 
of MW to be zeros we easily get the nearby rational positive semidefinite 
matrices M^'' for / = The nearby rational positive semidefinite 

matrices pl^.^l, Ql^d] 

can be recovered similarly. 
Having (^^(x) = v/-r^(x) and M H , . . . , Ql^'^l for 1 < / < g, . . . , 1 < j < p, 
the program is converted to 

inf Tr(MM,l^[«'oi,p[^'Ol,g[^'Ol)^ 

s.t. (^,„(x) - mKx)^ ■ MW ■ mKx)^Kx) 
^mo(x)^-MM ■mo(x), 

^e'{x.') - X;-=i w«/,i(x)^ ■ W't"''*] ■ wu',i{^)gu',ii^) 

- Yli=i v«',«(x)^ ■ ■ v«/,„(x)p«/,„(x, x') 

^w«,o(x)^-iy[«J°]-w«,o(x) 

T . p[£fl] 



P^.oixj 



q^,i(x)0,i(x) + e£,2 



^ -q,,o(x)^ ■ Ql^'Ol ■ q,,o(x) 

Observing ([HD, the matrices M'^^, . . . , Q^^'^^ have floating point entries, while 
the matrices Aft'l , . . . , Q^^'-'^ are rational positive semidefinite matrices. There- 
fore, the remaining task is to find nearby rational positive semidefinite ma- 
trices M^^\ . . . , Q^^'^^ such that the constraints in ([HD hold exactly. To fulfil 
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this task, we can first apply Gauss- Newton iteration to refine M'"', . . . , Q^^'^\ 
and then recover the rational positive definite matrices M^^\ . . . , Q^^'^'^ re- 
spectively from the refined M^^\ . . . , Q^^'^\ by orthogonal projection if the 
involved matrix is of full rank, or by rational vector recovery method other- 
wise. ^ 

Finally, we check if all the matrices M^, . . . , Q'^'"' are positive semidefi- 
nite. If so, then return (^^(x) > as an invariant of the given hybrid system 
at location £ G L; otherwise, return "we cannot find invariants of the given 
degree bound" . 

Remark 3. The above technique based on SOS relaxation and exact polyno- 
mial recovery can be applied to computing the inductive invariants of hybrid 
systems, which guarantee the safety of the given hybrid system. 

3.4- Algorithm 

The discussion in Section 13.31 leads to an algorithm of computing the 
(inductive) invariants of polynomial hybrid systems. As stated above, we 
only present how to compute the invariants ipi{'x.), for i E L, that satisfy (E]), 
and the case of computing the inductive invariants is similar. 

Algorithm Polynomial Inequality Invariant Generation 

Input: *■ H : (V, L, T, 6, V, io) a polynomial hybrid system. 

*■ d E Z>o: the degree bound of the candidate polynomial invari- 
ants. 

*■ D E Zi>o: the bound of the common denominator of the coeffi- 
cient vector of the polynomial invariants. 

*■ e E Z>o: the degree bound 2e of the SOSes used to construct 
the SDP system. 

*• T G M>o: the given tolerance. 
Output: *• <^^(x) > 0: the verified polynomia invariant at each location 
dEL. 

1. Compute the candidates of polynomial invariants 

(i) For each locaiton i E L, predetermine the templates of (/?^(x), with 
degree d, and construct an SDP system of the form ([H]), where the 
degree bounds of all the involved SOSes are 2e. 
• If the SDP system has no feasible solutions, 

return "we can't find polynomial invariants with degree < d 
at each location"; 
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• Otherwise, 

obtain a numerical vector v^, numerical constants e^^i, e£_2 and 
numerical positive semidefinite matrices M^'l, W^^^ yV^M^ 
QliJ] ioT < I < q,0 < i < s,l < u < t,0 < k < r,0 < 

J < P- 

(ii) For the common denominator bound D, compute from v| a ratio- 
nal vector \i by Diophantine approximation algorithm, and get 
the associated rational polynomial (^^(x). Similarly, the nearby 
positive contants and are obtained. 

(iii) Convert all the M^, . . . , Q'^'-'' into rational and positive semidefi- 
nite matrices M^, . . . , Q^^'^\ for 1 < / < g, . . . , 1 < j < p. 

2. Compute the exact SOS decomposition 

(i) Reconstruct an SDP system of the form to get approximate 
positive semidefinite matrices M^, . . . , Q'^'*^' satisfying ([H]). 

(ii) Apply Gauss-Newton iteration to refine the matrices M^^^ , • • • , 
obtained in Step [2] (i). 

(iii) From the refined M^^\ . . . , Q^^'^\ compute the rational matrices 
M^^\ . . . , Qt^'*^! respectively by orthogonal projection method if 
the involved matrix is of full rank, or by rational vector recovery 
if the matrix is singular. ^ ^ 

(iv) Check whether all the matrices M^"! , . . . , Q'^'"' are positive semidef- 
inite. 

• If so, return ^^(x) > as an invariant at location £ G L; 

• Otherwise, 

return "we can't find polynomial invariants with degree < d." 

Remark 4. Our algorithm cannot guarantee rational solutions will always 
be found since there exists limitations in the above algorithm on choosing 
the degree bound e and the common denominator bound D. Furthermore, 
it is difficult to determine in advance whether there exists invariants with 
rational coefficients or not. Therefore, even if our algorithm cannot find the 
invariants, it does not mean that the given hybrid system has no invariants 
with the given degree bound d. 

4. Experiments 

In this section, some examples are presented to illustrate our method for 
safety verification of hybrid systems. 
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Example 1. 



'2S 



, Example CLOCK] Consider a nonlinear continuous system 



X 




' -'iy + y'' 






Qx — x"^ 



with location invariant \l/ = {(x,?/)G]R^:l<a:<5Al<j/<5}. The 
problem is to verify that all trajectories of the system starting from the initial 
set B = {(x,?/) E M? : A < X < 4.5 /\y = 1} will never reach the unsafe 
set Xu = {{x,y) G : 1 < X < 2 A 2 < y < 3}. The safety of the 
continuous system can be verified if we can find a polynomial ip{x,y) which 
satisfies conditions in Theorem\^ We rewrite Q^X^-,^ as follows 

e = {(x, y)eM?: Oiix, y) > A O^^x, y) > A 9s{x, y) > 0}, 
^ = {ix,y) e : t/j,{x,y) > OA^x^y) > 0}, 
Xu = {(x, y)eR'': Ci(x, y)>OA C2(x, y)>0}, 

where 



6'i(x,y) = (4-x)(x-4.5), 62{x,y)=y-l, 93{x,y) = l-y, 

y) = {1- x){x - 5), ip2{x, y) = {1- y){y - 5), 
Ci(x, y) = (1 - x){x - 2), C2(x, = (2 - y)iy - 3). 

Assuming deg{Lp{x,y)) = d, for d = 1,2, ... and the degree bound of all the 
involved SOSes in the program ^ is 2e = 10. Then the SOS program ^ 
becomes 

ip{x,y) = aoix,y) + ai{x,y)9i{x,y) + a2(x, ?/)6'2(x, y) + as{x,y)9s{x,y), 
ip{x, y) = 0o(x, y) + 01 (x, ?/)V'i(x, y) + Mx, y)ip2{.x, y) + ei, 
-V9(x, y) = /io(x, y) + /ii(x, 2/)Ci(x, y) + /i2(x, ?/)C2(x, y) + £2, 

where ai{x,y),4)j{x,y), fik{x,y) G S2,2e,ei,e2 ^ We apply the algorithm 
in Section and increment d by 1 from 1 to 10 until a feasible solution of 
the SDP system is obtained. When d = 4, we obtain a feasible solution of 
the associated SDP system. Here we just list one approximate polynomial 

(^(x, y) = -4.3296 - 1.2975x - 0.10418y + 0.92562x2 + 0.18428xy 

+ 0.35738y2 + ... + 0.94032 x lO'^x^ + 0.17047 x 10"^?/^ 
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Let the tolerance r = 10~^, and the bound of the common denominator of 
the polynomial coefficients vector be 1000. By use of the rational SOS re- 
covery technique described in Section \3.3[ we obtain all the corresponding 
polynomials with rational coefficients, for instance, 



^, , 4113 1233 99 

ipix.v) = X y-\ X' 

950 950 950^ 950 



879 2 34 2 



7 



6 



95 



V 



38^^ 475 



xy 



46 
475^ 



:X 



Furthermore, a certificate of SOS representation shows ^p{x, y) satisfies the 
conditions in Theorem exactly. Therefore, the safety of this continuous 
system is proved. 



Example 2. 123 . Example ECO] Consider a predator-prey hybrid system 
depicted in Figure [H where 



/l(x) = /2(x) 



-Xi + XiX2 
X2 - X1X2 



g(l,2) : (X2 - 0.875) (12- 0.9) <() 




p(2, 1) : (.Ti - 0.7)'-' + (:i-2 - 0.7)2 < (, Qi 

Figure 2: Hybrid system of example 2 
The system starts in location £1, with an initial state in 

e = {{xi, X2) e : (xi - 0.8)2 _^ (^^^ _ Q_2)2 < q q^|_ 
Our task is to verify the system never reach the states in 

Xuih) = {{xi,X2) e : 0.8 < Xi < 0.9 A 0.8 < X2 < 0.9}. 

To verify the safety of this system, we need find the corresponding invari- 
ant polynomials ipi{xi,X2) and ip2{xi, X2) at locations £1 and £2, respectively. 
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Similar to Example 1, we construct the associated SOS system, and find 
the feasible numerical solutions from SDP solver: 

= 0.34871 - 0.45903x1 + 0.018001^2 + 0.2212x^ - 0.45764xiX2 + 0..17991x^, 
ip2 = 0.011167 + 1.2891x1 + 0.56568x2 + 0.88855x? - 0.56553xiX2 - 0..18386x^. 

Let the tolerance r = 10~^, and the bound of the common denominator 
of the polynomial coefficients vector be 1000. By use of the rational SOS 
recovery technique , we obtain all the corresponding polynomials with rational 
coefficients. The invariant polynomials with rational coefficients are 



(/?l(xi,X2) 



329 433 17 

^1 + 777T^2 



27 



209 



944 944 
11 1217 



944 944 



■Xi 



944 
267 

472 



-X1X2 + — Xi + 

267 



85 
472' 



-X 



2i 



X2 



472 



839 2 87 2 



Furthermore, all the remaining related polynomials in ^ can be written as 
SOSes of the polynomials, which means tpi and ^2 satisfy all the conditions 
in Theorem [31 So the safety of hybrid system is proved. 



s(l,2):(j-2-1.6)(.-r2-2)<0 




p(2,l) : (.1-1 - 1)2 + (.T., - \f < 0.01 



Figure 3: Hybrid system of example 3 



Example 3. Consider a hybrid system depicted in Figure where 

2xi + X2 



/i(x) 



/2(X) 



X2 

-Xi + X2 



X1X2 — X2 — 1 

The system starts in location ii, with an initial state in 

e = {(xi, X2) e : (xi - 1.5)2 + xl< 0.25}. 
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Our task is to verify the system never reach the states in 



Xu{ii) = {{xi, X2) e W : (xi + 1)" + (X2 + ly < 0.16}. 

To prove the safety of the hybrid system, it suffices to find an inductive 
invariant polynomial ip{xi, X2) which satisfies all the conditions in Theorem\^ 

Using the same techniques illustrated in Examples 1 and 2, we obtain the 
inductive invariant polynomial with rational coefficients 

22 319 251 239 2 

(piXi, Xo) = 1 Xi Xo H X, . 

' 49 931 931 931 ^ 

Moreover, satisfies the conditions in Theorem\^ exactly. Therefore, the in- 
ductive invariant can guarantee the safety of the hybrid system. More details 
about the verification of conditions in Theorem\^ based on SOS representa- 
tions of polynomials with rational coefficients can be found in Appendix. 



5. Conclusions 

In this paper, we present a symbolic-numeric approach to compute in- 
equality invariants for safety verification of hybrid systems. Employing SOS 
relaxation and rational vector recovery techniques, it can be guaranteed that 
an exact invariant, rather than a numerical one, can be obtained efficiently 
and practically. This approach avoids both the weakness of numerical ap- 
proaches to verify safety of hybrid systems and the high complexity of sym- 
bolic invariant generation methods based on quantifier elimination. 
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6. Appendix 

Solution to Example 3: 

The initial state 9, the unsafe region Xu{ii), the state invariant "^{i), 
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the guard condition g{i,i') and the reset map p{i,i') can be expressed as 



e = {{x^,X2) G : ei{xi,X2) > 0}, 

^(1) = {(xi, X2) e : tpi,iixi, X2) > A ^i,2(xi, X2) > 0}, 
*(2) = {{xuX2) e : iJ2Axi,X2) > A ^^2,2(3^1, Xs) > 0}, 
(7(1,2) = {(xi,X2) G : ^i2(xi,X2) > 0}, 
(7(2, 1) = {(xi, X2) G M^' : g2i{xu X2) > 0}, 
p(l,2) = {(xi,X2) G : Pi2(xi,X2) > 0}, 
p(2, 1) = {(xi, X2) G : P2i(xi, X2) > 0}, 
X„(l)(xi,X2) = {(xi,a;2) G : Ci(a;i,a;2) > 0}, 

where 

e^{xi,X2) = 0.25 - (xi - 1.5)2 - xl, 

Ci(xi,X2) = 0.16 - {x^ + 1)2 - (x2 + 1)2, 

^1,1(2:1, 2:2) = (2:1 + 1)(2 - Xl), ^i,2(a;i,X2) = (xg + 1)(2 - X2), 

V'2,i(2;i,a;2) = (xi - 2.5)(3 - xi), ^/'2,2(xi, X2) = (x2 - 2.5)(3 - X2), 

^12(3:1, X2) = (x2 - 1.6)(2 - X2), fi'2i(xi,X2) = (x2 - 2.5)(2.75 - X2), 

Pi2(a;i, X2) = 0.01 - (xi - 2.6)2 - (x2 - 2.8)2, 

P2i(a;i,X2) = 0.01 - (xi - 1)2 - (x2 - 1)2. 

Let the tolerance r = 10~2^ and the bound of the common denominator 
of the polynomial coefficients vector be 1000, we can find that the inductive 
invariant polynomial ^(xi,X2) satisfies 

^(xi,X2) = ao(xi,X2) + ai(xi,X2)6'i(xi,X2), 

^(Xi,X2) = Ai2o(2:i,X2) + Ai2l(Xi,X2)^l2(Xi,X2) + 7l2 (^^l , X2)pi2 (a^l , 2:2) , 

^(xi,X2) = A2io(2:i,X2) + A2ii(xi,X2)fi'2i(a;i,X2) + 721(3^1, a;2)p2i (2:1, 2:2), 

^(Xi,X2) = 01o(2;i,X2) + 01l(xi,X2)V^ll(Xi,X2) + 012 (2:1 , X2)^12 (Xl , X2) + €1, 
^(Xi,X2) = 020(2:1, 2:2) + 02l(2;i,X2)V^2l(2;i,X2) + 022 (2:1 , X2)l/'22 (2^1 , X2) + £2, 
-^(Xi,X2) = Po(Xi,X2) + /il(Xi,X2)Cl(2:i,X2) +6, 
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where 



i r)( rp rp \ ^2 I 319 rp ^51 ^ _|_ 239 ^2 

V^V-^l) •>^2} 49 931 1 931 2 "T 93]^ -i' 



^y^X-, I] 49 ^ 931 1 931 2 ~ 931 1' 

n — 838 _ 1565™ _ 251^ , 867^2 , 628^2 ~ _ 628 

uol,xi,X2; 93;^ 931-^1 931-^2 -r 931-^1 -r 931-^2 > ^1 931 

X (t T \ — — 4463^ _ 14169 ^ , 472 2 , 12 ™2 

Ai20l,'i'l, -^2; — 1900 4655-^1 4655 "^2 ^ 931-^1 19-^21 

A — 355 ~ _ 233 

^121 931 ) 712 931 ) 

\ irp rr \ — 133187 _ 1271„ _ 1997^ , 1034^2 , 1110^2 

/V2l0l,-Al,"L2j — 37240 931-^1 532-^2-1- 931 -^1 "T 931-^2) 

^211 = 1331 721 — 931) 

X ^ ^ — 153 _L I 298^ I 227^ ^ , 971„2 , 272 2 

t'IOV'^I ) '^2J 931 ~r 133-^1 ~r 931 ■^2 ~r 93i'^l-t'2 ~r 93i'''l ~r 931 ■''2' 



h = ^ rh = ^ 7 = 26 

Illl 931' ^12 133, tl gg-^, 

J, ,r \ — 3751 _ 3337^ _ 3373™ , 349 ^2 , 478™ ™ , 319 ™2 

^201-^1 ) -^2^ 931 1862 1 1862 2 "T 931 ■''1 ~r 93i"''l-''2 "T 931-^2' 

'^21 931 ) 931 , C2 931 , 



/io(a;i, X2) — ~ 97-^1 ~^ "97^ '^2 + QjXlX2 + gy^^i + 

~ 564 ~ 58 

"1 ~ 931 ' 931 ■ 



97-^2' 



The exact SOS representations of above polynomials are as follows: 

-X (rr ^ \ — 931 l2 I 3120712 ^2 , 380230641 ^2 
Oo{-l-i,X2) — 838"'ll 2042055'''12 46469677 "■13' 

T /'^ ™ ^ _ 1900 ^2 I 127778819 l2 , 12831251475 l 2 
^120l-^li-^2j — 8403'''21 "T 13782225 "'22 "T 2601209876 "23' 

T 1'^ ™ ^ _ 37240 l2 _|_ 991976776 l2 , 38289861701 ^2 
^2101-^1, -^2; — i33i87"31 205638355 "32 13835779654 "33' 

I f™ „A _ 931/,2 I 142443 ;^ 2 , 72301460 ^2 
</^10l-^, — 153 "-41 19415 "42 4096293 "43' 

J _ 931 ^2 I 55874896 1,2 , 7231984725 ^2 

<P20l'i, - 3751 ''-51 + 7767975 "52 + 1110827906 "53' 

77 1'™ — 3325 7,2 I 978432 u2 , 29133446909 1,2 

/i-Ol-^1, -^2; — 4992"'61 235283"62 941925575 "63' 



where 



h = 838 _ 251 _ 1565 L _ 2042055 _ 392815 L _ 46469677 

"11 931 1862 2 1862*^1' "12 3120712"^2 3l20712"''l' "13 380230641 '^l ' 

L _ 8403 _ 14169 _ 4463 _ 13782225 _ 21078749 L _ 2601209876 

21 1900 9310 2 9310-^15 "•22 127778819 2 127778819 1' 23 12831251475 1" 

u _ 133187 _ 1997 _ 1271 _ 205638355 _ 12690935 L _ 13835779654 

31 37240 1064 2 1862 1' 32 991976776 2 35427742 1' 33 38289861701 1' 

h = 153 _L 149^^ _L h = 19415 _ 29048 l ^ 4096293 

"41 931 ~'~ 931 2 "T 266 1 ' 42 142443 2 142443 1' 43 72301460 1' 

u _ 3751 _ 3373 _ 3337 l _ 7767975 , 3088123 l _ 1110827906 

"51 931 3724-''2 3724-^1' 52 55874896 2 "r 55874896 1 ' 53 7231984725 1' 

1^ _ 4992 I 197 I 809 l _ 235283 _ 3984325 l _ 941925575 

"61 3325 ~'~ 266 2 "r 1862"''1' 62 978432"^2 18590208 1' 63 29133446909 1 ' 
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